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Abstract 

We study the use of effective transfer matrices for the numerical computation of 
masses (or correlation lengths) in lattice spin models. The effective transfer matrix 
has a strongly reduced number of components. Its definition is motivated by a 
renormalization group transformation of the full model onto a 1-dimensional spin 
model. The matrix elements of the effective transfer matrix can be determined by 
Monte Carlo simulation. We show that the mass gap can be recovered exactly from 
the spectrum of the effective transfer matrix. As a first step towards application 
we performed a Monte Carlo study for the 2-dimensional Ising model. For the 
simulations in the broken phase we employed a multimagnetical demon algorithm. 
The results for the tunnelling correlation length are particularly encouraging. 



CERN-TH.6844/93 
March 1993 



1 Introduction 



The computation of the mass spectrum is one of the major goals of modern quantum 
field theory. In the Euclidean path integral formulation on the lattice PJ, the mass 
spectrum can be recovered from the large-distance behaviour of suitably chosen 
correlation functions. The correlation functions can (in principle) be computed 
with the Monte Carlo method ||. 

The bridge between the path integral and the Hamiltonian formulation is in the 
transfer matrix [ED] . The mass spectrum can be directly read off from the eigenvalues 
of the transfer matrix. The task is therefore to diagonalize the transfer matrix. An 
exact solution of the full problem is possible only in a restricted class of models. 
The most prominent example for a model that can be solved via the transfer matrix 
approach is the 2-dimensional Ising model ||. A direct numerical diagonalization 
of the transfer matrix is restricted to systems with very small spatial extension [|]]. 

As the origin of the problem with the transfer matrix approach is in the huge size 
of the matrix to diagonalize, it is natural to reduce the number of degrees of freedom 
by some 'coarse- graining procedure'. Inspired by block spin renormalization group 
ideas ||, we define an 'effective transfer matrix' .[] This effective transfer matrix 
has a drastically reduced size and can be diagonalized with standard numerical 
procedures. 

In section || we motivate the concept of the effective transfer matrix, starting 
from a block spin renormalization group point of view. We arrive at an intuitive 
derivation of three different 'rules' for a definition of an effective transfer matrix. In 
section |^ we show (for one of the rules) that the mass of the first excited state can 
be exactly recovered from the effective transfer matrix. To learn more about the 
properties of the effective transfer matrix, we then study the 2-dimensional Ising 
model on small lattices. Here we can compute the spectrum of the effective transfer 
matrix without any approximation and compare it with the exactly known spectrum 
of the full transfer matrix. So we can study to what extent the low-lying spectrum 
of the full transfer matrix can be recovered from the effective transfer matrix. In 
the final section we present results for the low-lying spectrum of the 2-dimensional 
Ising model on lattices up to 64 x 2048 as obtained from effective transfer matrices 
computed by Monte Carlo simulations. We discuss results from simulations at the 
critical point and in the broken phase of the model. The multimagnetical demon 
algorithm we employed for the simulations in the broken phase provided us with high 
statistics for tunnelling events. The success is based on an increased probability 
for configurations with low magnetization and a high performance due to a very 
efficient implementation. The longest tunnelling correlation length we measured 
was 43500 ± 1600 on a 64 x 128 lattice, consistent with the exact value 44014.4 



1 In the literature, the notion of an 'effective transfer matrix' is used in a variety of different 
contexts. The reader might be interested to compare our notion with that developed, e.g. in Q, 
|?J and ||. Furthermore, for one of our definitions of the effective transfer matrix, our approach 
is very similar to a method that has been discussed in the context of glueball mass calculations 
||. [T(|. See also Appendix C. 



2 Effective 1-Dimensional Model 



We consider a (d+ l)-dimensional Euclidean quantum field theory on a cubic lattice 
A of size N% x N T . Here, N a denotes the spatial extension of the lattice, and N T 
is the extension in the time direction. Let us label the sites of the lattice by (t,x), 
where t runs from 1 to iV T , and the spatial coordinates Xi cover the range 1 . . . N a . In 
the Euclidean path integral formulation one integrates over stochastic variables <pt,x 
that are attached to the sites of A. The model is defined by the partition function 

Z = £exp(-W(0)). (1) 

Let us denote the configurations of on 'time slices' t by <p t . Assume that H(<f) 
takes the following form: 

= + (2) 

t 

Then for periodic boundary conditions in the t-direction the partition function may 
be written as 

Z = Tr(T N r). (3) 
The transfer matrix T can be chosen to be symmetric: 

T(<p, <p') = exp(-K(<p, <pf) - \[V{y) + V(<p')}) (4) 

<f and ip' are again time slice configurations, i.e. configurations on rf-dimensional 
sub-lattices of size N%. We shall also use a 'bra' and 'ket' notation, e.g. 

r( ¥> ,^) = ^|r| v />. (5) 

Let Aj, i = 0, 1, . . . denote the eigenvalues of T, such that A > Ai > A 2 > 

The corresponding eigenvectors shall be denoted by For eigenstates with zero 
momentum, masses are defined as 

mi = - In ^ = \ . (6) 

^0 Si 

The £j denote the corresponding correlation lengths. We shall now consider block 
spin transformations that transform the statistical system from a (d+ l)-dimensional 
one to a 1-dimensional one. Block spins are defined as averages over time slices: 

$ t = 7V-^^. (7) 



X 



For integer I > 1 we define an 'effective Hamiltonian' 

exp(-w2(<&)) = ^exp(-W)) I] K*t ~ ^EM » ( 8 ) 
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with Gi = {1, 1 + 1, 21 + 1, . . . , N T — I + 1}, and Z can be rewritten as 

z= (nz) «p(-«2(*)) . (9) 

Let us now relabel the sites in Gi with integers that run from 1 to N' T = N T /l. Let 
us furthermore assume that (after the relabelling) for some / with good precision 
7"^g($) couples only nearest neighbours, i.e. is of the form 

w2(*) = E (*2(*«, $ m) + • (io) 

Note that, in general, the effective Hamiltonian will contain interaction terms of 
arbitrary range. However, one can generally assume that these terms decay expo- 
nentially with the distance and might be neglected for a first qualitative analysis. 

Assuming a nearest neighbour interaction effective Hamiltonian, we can again 
rewrite the partition function in terms of a transfer matrix 

Z = Tr([T^), (11) 

and 

TjJ($, $') = exp $') - \\V$(*) + V® (<&')]) ■ (12) 

What have we gained? The number of degrees of freedom has been considerably 
reduced. Take for example the Ising model. Here the original transfer matrix is 
2 N ° by 2 N °. If we choose as block spin the time slice magnetization then is a 
matrix with + 1 by N% + 1 components. So there is a drastic simplification in the 
eigenvalue problem. The crucial question is, of course, to which extent the spectrum 
of T e fj reflects properties of the original system. This question will be studied in 
detail in the next section. 

Let us first discuss how one can express the elements of as expectation 
values in the statistical system. This is necessary if we want to compute them, 
e.g. by Monte Carlo simulation. To simplify the discussion, we shall assume in the 
following that the $ t take discrete values only. Let us denote the discrete values of 
$t by M,N, .... They correspond to 'improper states' \M)),\N)), etc. We shall use 
throughout 'double brackets' to denote states in the space that Tjg acts on. We 
now introduce the operator 5$ tt M that takes the value 1 if $ t = M, and 5$ 4) m = 
elsewhere. Assuming periodic boundary conditions in the t-direction, the correlator 
of two such operators at distance one can be written as follows: 

In the following we shall sometimes omit the extra factor Z -1 , because it leads only 
to an irrelevant shift of the ground-state energy. Masses and correlation lengths are 
unaffected. Equation flT3D has to be resolved with respect to the effective transfer 
matrix elements. This is particularly easy for the case N' T = 2, 

((M\T®\N)) = y/<5* uM 5* 2 ,N>. (14) 
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We say that the effective transfer matrix is defined according to the symmetric 
periodic lattice rule. For larger N' T (asymmetric periodic lattice rule), the solution 
for can be found by iteration (cf. Appendix A). The situation becomes simple 
again in the limit N' T — > oo, 

< S* t , M 5* t+1 , N >= ((0\N))((N\TS\M))((M\0)) . (15) 
Using the fact that < 5$ uM >= ((0|M))((M|0)), one finds 

mT m m = <Ww>> (16) 

\J < 0<s> u M >< 0<s> u N > 

We say that is defined here via the N T = oo rule. 



3 How the Mass Gap is Recovered 

In this section we want to discuss how the effective transfer matrix is related to 
the transfer matrix of the basic system, T, and how much of the spectrum of this 
matrix can be recovered from the spectrum of the effective transfer matrix. We will 
concentrate on the discussion of the N T = oo rule. 

We define 'magnetization pieces of the ground-state' by 

\M)=U M 1 Y.^V = M) 1^10). (17) 

f 

Here, \ip) are time slice configuration states, and 

y? = iV- d X>. (18) 

X 

The constants Mm are chosen such that (M\M) = 1. One can convince oneself that 
(N\M) = 5n,m- However, the |M)'s do not, of course, span the complete Hilbert 
space. Let us denote by V the projector onto the subspace spanned by the |M)'s , 

V = J2\M)(M\. (19) 

M 

The effective transfer matrix for the N T — oo rule is given by 

iM\T®\N)) = (M\T l \N) . (20) 

Let us denote the eigenstates of by and the corresponding eigenvalues by 
A?. We now claim that A^ = Ag. The proof is as follows: 

<<M|T«|0>> = £<< M l T cS|A0><<^|0>> 

N 

= y,{ m \ t1 \ n )( n \°) 

N 

= (M\T l V\0) 

= A ((M|0)). (21) 
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Here we have used the fact that V\0) = |0). 

What about the first excited state? Let us try the following wave function: 

({M\V))=c- 1 ll2 {M\l). (22) 

The constant c\ is chosen such that this state is normalized to 1: 

Ci = E( M l 1 ) 2 = ( 1 ^ 1 )- (23) 

M 

The action of the effective transfer matrix on this state is easily found to be 

((M\T$\1')) = Ci 1/2 (M\TV\1) . (24) 

Now consider the first excited state 1 1) of the original transfer matrix T with eigen- 
value Ai and rewrite it as follows: 

\l) = <h 1 P\l)+ y £fi,M- ( 25 ) 

i>l 

That only states with i > 1 contribute in the sum is due to the fact that (0 17^| 0) = 1, 
and (0[7- > | 1) = 0. Solving this equation for V\l) and inserting this in eq. (|24|) , one 
obtains 

{(M\TS\1>)) = Cl \[ ({{M\l')) - E^ (^) W>) • (26) 
For large / one therefore has 

TiSn = c 1 \[\n+o^y^ . (27) 

From the spectrum of the effective transfer matrix, we obtain estimates for the mass 
of the first excited state 

m? -HSH- (28) 

It is easy to see that mf behaves like 

mf = mi y~ + exponentially small corrections . (29) 

One can get rid of the 1/7 corrections by combining pairs of results with different 
Vs. For the correlation length of the first excited state one gets 

£i = — — h exponentially small corrections . (30) 

h — h 

What about the higher states? Let us first see what happens with the second 
excitation. Proceeding naively, we would start with 

<<M|2'>) = C - 1/2 <M|2>, (31) 

and go along the same lines as in the proof above. We would then find that every- 
thing would work just the same as for the first excited state, provided that (1|"P|2) 
vanishes. (This can happen for symmetry reasons, see the discussion of the Ising 
model below.) Generally, when for higher states the overlaps (i\V\j) are very small 
for i<j, then the state allows the reconstruction of Aj with good precision. 
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4 Exact Results on Small Lattices 



As a first study of the effective transfer matrices defined above, we made numer- 
ical calculations for the 2-dimensional Ising model on lattices with small spatial 
extension. These calculations are exact in the sense that only standard numerical 
techniques were used. No Monte Carlo simulations were involved. 
The model is defined through its partition function 

Z = exp /3 a * a y ■ ( 32 ) 

cr x =±l \ (x,y) J 

Here we have denoted lattice sites by x = (t, i), and (xy) denotes a nearest neighbour 
pair. The (symmetrically chosen) transfer matrix of this model is given by 

T(<f> t , (fit+i) = exp - ^ Vt,i°t,i+i + P X] a t,i<7t+i,i + crt+i,t^t+i,i+i ! ( 33 ) 

i=l i=l ^ i=l ) 

4> here denotes time slice configurations of the Ising spins. We assume periodic 
boundary conditions in the space direction, i.e. the site with i = N a + 1 is identified 
with the site i = 1. 

The 2-dimensional Ising model was first solved by Onsager in 1944 ||. The 
eigenvalues of the transfer matrix can be read off easily from ||11|| . See also the 



reviews in ref. 12 



As a consequence of the fact that the transfer matrix commutes with the spa- 
tial translation operator, the eigenstates can be chosen as simultaneous eigenstates 
of energy and momentum. They can further be classified in states that are sym- 
metric/antisymmetric with respect to reversal of all spins. We shall denote the 
corresponding masses or correlation lengths with a subscript s or a, respectively. 

In the infinite volume limit, i.e. when N a — > oo, the spectrum is twofold de- 
generate in the broken phase {(3 > (3 C ). For finite N a the degeneracy in the low- 
temperature regime is lifted and level splitting occurs because of tunnelling. 

We will consider the time slice magnetization as effective spin in the following. 
In this case only the zero-momentum eigenvalues can be recovered from the effective 
transfer matrix. f\ 

We shall give results for the N T = oo rule and for the symmetric periodic lattice 
rule. We have results for N a < 8. The basis for our computations is an accurate 
determination of the eigenvalues and eigenfunctions of the transfer matrix defined 
in eq. (|33|). In principle the eigenfunctions of T are given in the work by Kaufmann 



1341 - However, the expressions seem too complicated to be of practical use for our 
goals, so we decided to use the computer. 

Let us first discuss the case N T = oo. We numerically determined the eigen- 
states \i) and the 'magnetization pieces' \M) defined in eq. QT7p. Together with the 
eigenvalues of T, they enter the matrix elements of Tjjj as follows: 

((M\T^\N)) = Y,A(Mm\N) ■ (34) 



2 Note that effective states that couple to nonzero momentum could also be employed 
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The |M)'s are zero momentum states. Hence we only needed to calculate the trans- 
lational invariant eigenstates of T. 

In the case of the periodic system with N T = 21 (symmetric periodic lattice rule) 
one finds 

({M\T®\N)) = ^^(i\M)(M\j){j\N)(N\i). (35) 

For both rules we computed the spectrum of T e g for various 

/ and determined the £™ defined in eq. fl28[) for the lowest-lying states. Our 
results for the N T = oo rule are listed in tables [1| and The results for the symmetric 
periodic lattice rule are quoted in tables |] and £|. We always quote the results for 
I — 1, 2, 4 and 8. In addition, we give an estimate for the 'true' £ obtained by 
combining ^ and ^ according to eq. (|3~0"D. The exact results are always quoted 
with an 'e' in the second column. 

Let us first look at the results for the N T = oo rule. As we proved in section |3|, the 
£o,a should converge towards the exact limit £ ,a- This can indeed be observed. The 
second largest correlation length £ 1)S also nicely converges. This is a consequence 
of the fact that the corresponding state is symmetric in the magnetization while 
the first excited state is antisymmetric. Therefore the matrix element (lsl^lOa) 
vanishes, and the proof that the correct correlation length can be recovered goes the 
same way as for the first excited state. The convergence for £i ja is also quite good. 
This is anticipated from the discussion at the end of section |3|, since c 0a is very close 
to 1. Following the discussion of section the constant q = Sm(^K) 2 generalizing 
eq. (|23|) should indicate the goodness of the effective eigenstate corresponding to 
Ci = 1 indicates a perfect representation of the eigenfunction \i) by the corresponding 
effective eigenfunction 

For L = 4, 6, 8 and 10 we calculated this quantity explicitly, starting from the 
exact eigenfunctions with labels i = (0, a), (l,s) and (l,a). The results for the 
/?-range 0.35 - 0.55 are given in fig. 1. We make the following observations: The 
overlap co a is very close to 1. It becomes better for (3 — > and for (3 — > oo. The 
overlaps c\ s and c\ a are worse. The deviation from 1 is in the range of several 
per cent; C\ s and c\ a are 'best' in the critical domain. All c's become worse with 
increasing N a . It is an open question whether they have a finite limit in the N a — > oo 
limit. 

The t^ps do not converge. Here the states of the original system seem to mix 
somehow to the effective states. 

Let us now turn to the discussion of the results for the periodic symmetric rule. 
Table ^| is encouraging: when we combine the I = 4 and I = 8 data, all correlation 
lengths are reproduced within the given numerical accuracy, with the exception of 
^2,s- However, this is only so on very small lattices. A careful look at table (| 
reveals that only the leading correlation lengths £q a converge to the exact value 
£o.a- The estimates for the other £'s are not entirely off. They might be regarded as 
unsystematic approximations. However, we do not know what happens if the spatial 
extension of the lattice is further increased. See also the discussion of Monte Carlo 
results on larger lattices in section [|. 
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5 Monte Carlo Results for the 2D Ising Model 



5.1 The Critical Model, N T = oo rule 

At the critical coupling (3 C = 0.4406868... of the 2-dimensional Ising model, we 
performed Monte Carlo simulations on N a x N T = 16 x 512, 32 x 1024 and 64 x 2048 
lattices. We always measured and stored the magnetizations of all time slices after 
five Swendsen-Wang updates ]13|] of the entire lattice. 



We made 5000 measurements for N a = 16, 3900 measurements for N a = 32, and 
3000 measurements for N a = 64. 

For the simulations at criticality we used the N T = oo rule to get estimates for 
the effective transfer matrix elements. 'Infinitely long' lattice here means, of course, 
that N T is much larger than the largest correlation length involved, so that the 
influence of a finite lattice extension in the t-direction can be neglected. In all three 
cases, N T /^ 0ta is about 25. This is certainly on the safe side. 

The eigenvalues and eigenstates of the effective transfer matrix were then cal- 
culated using standard numerical procedures. From the eigenvalues of the effective 
transfer matrix we determined the for various distances 

In order to estimate the statistical errors of our results we divided the whole 
sample in 5, 10, 20 and 40 bins and computed the effective transfer matrix and its 
eigenvalues separately within each of the bins. The statistical errors quoted in the 
tables were obtained from the mean square fluctuations over the outcomes of the 
various bins. We considered the error estimate as reliable when it was approximately 
independent of the number of bins. 

This binning analysis was also used to check for a bias due to a too small statis- 
tics. We averaged the results from the bins and compared the outcome of various bin 
sizes. Only those results can be trusted where there is an agreement of the results 
over the various bin sizes. The quantities for which we obtained stable estimates 
are quoted in tables |], || and 0. 

With increasing I, at least the largest two correlation lengths £o,a and £i jS converge 
towards the corresponding exact value. Similar to the exact results on small lattices 
we find the estimates for £ l a converging within the statistical accuracy towards the 
exact results. 

The tables also contain our final estimates for the correlation lengths obtained 
by combining pairs of effective correlation lengths with different I according to eq. 
(PUp. The convergence and stability of these combined estimates shows that the 
I dependence is indeed of the form of eq. fl2"9"|). 

In fig. 2 we present plots of the effective wave functions (eigenfunctions of Tjj) 
for N ff = 64 and for / = 5, 10 and 20. One can observe that the effective wave 
functions obtained from different I are almost identical. For I = 20 the effective 
wave function for | la) becomes noisy because of insufficient statistics. 

In order to compare the effective transfer matrix approach with a 'traditional' 
determination of masses, we present here the results of a conventional evaluation 
of our Monte Carlo data for the 2-dimensional Ising model at criticality. As an 
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example we choose N T = 64 and N T = 2048. We define 

G 0>a {t) = < (p (pt > , 

G 1>s (t) = < <&<g > ~ < 0o >< <fi > • (36) 

These correlation functions behave like Gi(t) = g^ exp(— + ... . We define 
'effective correlation lengths' as 

ef' t = l/(ln(G i (t))-ln(G,(t + l))) . (37) 

With increasing t these quantities should converge to the true correlation lengths. 
Table || shows our results for the ^ S,t , for i = (0, a) and i = (1, s). We used exactly 
the same Monte Carlo data as for the effective transfer matrix. The comparison of 
the results obtained with the two different methods shows that the transfer matrix 
results have statistical errors that are roughly a factor of two smaller than the 
'conventional' ones. 

We evaluated the same Monte Carlo data with a third technique: Based on the 
results for the eigenfunctions of the effective transfer matrix T e ^, we determined 
observables that are expected to have improved overlap with the eigenstates |0,a) 
and 1 1, s) , respectively. For details see Appendix C. With the improved observables, 
we again computed the 'effective correlation lengths' ^ S,t , for % = (0, a) and % = 
(l,s). The comparison of these quantities, 'standard' and 'improved' is shown in 
figs. 3 and 4. At least for £ , a > the improvement is striking. 



5.2 Monte Carlo Results for the Broken Phase 

We also did simulations in the broken phase, at (3 = 0.47. This value is low enough 
so that the tunnelling correlation length becomes very large even with modest N a , 
but is close enough to (3 C so that the bulk correlation length (4.349 . . . when N a = oo) 
is still substantially larger than one. 

In order to fight the supercritical slowing down due to exponentially suppressed 
tunnelling rates, we employed a multimagnetical algorithm. Thanks to multispin 
coding and the usage of demons to implement the nonlocal magnetization-dependent 
interaction terms in the multimagnetical ensemble, the algorithm performed very 
well on a CRAY X-MP. For more details, see Appendix B. We performed simulations 
for lattices with spatial extensions N a = 16, 32 and 64. For an overview of the runs 
and the statistics see table |^. 

For the analysis in the broken phase, we used the asymmetric periodic lattice 



rule throughout. Our results are quoted in tables [10], |TT] and [12]. We performed the 
simulations with N T = N a /2, N a and 2N a . Generally, we observe two trends: the 
larger N T is, or the closer one is to the symmetric periodic lattice rule, the smaller 
the systematic error in the tunnelling correlation length measurements becomes. 
When N T = 2N a , we find an impressive reproduction of the tunnelling correlation 
length that gets as big as 44014 on the N a = 64 lattice. Note that the tables also 
include results for the symmetric periodic lattice rule: the last row in each block of 
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the table quotes the results for I = N T /2. We also present the tunnelling correlation 
length data for the N a = 64 lattices in fig. 5. 

We also tried to reproduce the large tunnelling correlation length on the N a = 64 
lattice using the standard technique of fitting the correlation functions of time slice 
magnetizations. However, this did not lead to any sensible result. To us this seems 
to be a situation where our effective transfer matrix technique is completely superior. 
We therefore consider our method well suited for the study of the interface tension 

El- 

The other correlation lengths besides the tunnelling length are only approxi- 
mately correct. It seems that they do not converge. Here certainly further investi- 
gation is needed. 

The effective wave functions for (3 = 0.47 and N a = 64 are plotted in fig. 6, 
for / = 8, 16, and 32. That the tunnelling length is so large corresponds to the fact 
that the square of the ground-state wave function and the square of the first excited 
state differ only in a small neighbourhood of zero magnetization. 

6 Summary and Conclusion 

We have studied the use of effective transfer matrices for the computation of masses 
of a Euclidean quantum field theory, or, equivalently, the correlation lengths of a 
statistical mechanical model. Many questions are left open and deserve further 
theoretical and numerical investigation. The theoretical study of the other rules 
besides the N T = oo rule is not yet complete. Furthermore, one could consider the 
usage of blocks with finite extension in the time direction. In principle, one can 
employ any sort of 'effective spin' (not just the magnetization or its absolute value). 
Application of the method to systems with continuous degrees of freedom would 
also be interesting. 

A report on a study of the 3-dimensional Ising model will be published elsewhere 

m. 
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Appendix A: Matrix Equation for Asymmetric Rule 

The problem of solving the equation 

< 5* t , M 5* t+uN >= ((M\tS\N})((N\[tS} n ^\M)} (38) 
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with respect to the components of the effective transfer matrix (given the left-hand 
side, e.g. as an outcome of a Monte Carlo simulation), may be stated as follows: For 
a given symmetric matrix B with matrix elements By, find a matrix A with matrix 
elements A^ such that the following equation holds: 

By = An {A n )ij , (39) 

where A n denotes the n th matrix power of A. We do not know of a closed solution of 
this problem, but, inspired by the famous iteration prescription for the problem of 
finding the square root of a number a, x — > \{x + -), we try the following iteration 
prescription for the solution of eq. (|39|) : 



^-T^T^ + CtSB; (4°) 



1 + C J (A n ] 



'.J 



where ( denotes a parameter that can be tuned in order to optimize the convergence. 
With the exception of a few special cases the algorithm converged, although the 
convergence rates were sometimes very slow. 



Appendix B: Multimagnetical Demon Algorithm 

We here present a brief description of the multimagnetical demon algorithm used 
for the simulations at /3 = 0.47. 

When (3 > /3 C , the tunnelling rate between the ordered states becomes exponen- 
tially smaller as the size of the system is increased. However, in order to extract 
the tunnelling correlation length we need a good statistics of the tunnelling events - 
i.e. the configurations with the magnetization between the bulk expectation values 



±AAb- The multimagnetical method |L5| solves this problem by artificially enhanc- 
ing the probability of these states. In the standard (no demons) approach, this 
is achieved by modifying the probability of the spin configuration a with an extra 
weight function G(M): 

P mm (a)^e^ n ^G(M,), (41) 

where is the magnetization of the configuration a. The probability of the 
magnetization M. becomes 

Pmm(M) oc Pcan(-M) G(M) oc f ^ e~ m ^5(M, - M)j G(M) . (42) 

We use a notation where the /3-dependence of G and P is suppressed. Usually one 
aims at a constant probability distribution between implying that G(Ai) cx 

l/-f'can(-A / l). However, since P C an is unknown, one has to use an approximate form 
instead, which can be obtained, for example, by scaling up the function G used in 
the simulations performed with smaller volumes. This is often further refined by 
performing test runs and adjusting G(M) until satisfactory Pmm(-M) is obtained. 
For further details, see refs. [TJ|, |T5| ]. 
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The function G(Ai) depends on the global magnetization, rendering the update 
non-local and thus preventing straightforward vectorization and parallelization. In 
this work modified the standard multimagnetical algorithm by utilizing demons to 
change the magnetization; this approach enables us to use powerful multispin coding 
and is highly vectorizable. This method is a magnetic analogue of the multicanon- 



ical demon algorithm presented in [fLq| . The actual update becomes a two-stage 
process: first, the demon magnetization K\ D is changed by coupling the demons 
to a "multimagnetical heat bath", and second, the spin system magnetization Ai,s 
is changed by coupling the spins to demons and performing the spin update while 
preserving the total magnetization M. T = Ai D + Ais- More precisely, we perform 
the simulation according to the joint probability distribution 

P(M S , M D ) cc P can (M s ) n D (M D ) W(M T ) , (43) 

where n^M.^) is the number of demon states with magnetization Md, and W is 
a new weight function that depends only on M. T . The canonical expectation value 
of an observable O can now be obtained by reweighting: 

< >= S,ew-'(Afc) 

where (9; and A4t,i refer to the individual measurements of the corresponding quan- 
tities. Summing over M. D in eq. ([43]) and comparing it with eq. (f42"D , we note that 
G and W are related: if G(A4s) J2m d n D(M-D) W(M.t)i we obtain a similar 
probability distribution for .Ms in the two cases. In this work we are using demons 
carrying ±1 units of magnetization; with Np demons, the demon density of states 
becomes 

N D \ 

Ud{Md) = [(N D + M D )/2MN D -M D )/2]\ ■ (45) 
Let us look closer at the individual update steps. 

• In the multimagnetical demon refresh step the new Ain is chosen with prob- 
ability 

P{M D ) (xn D {M D )W{M s + M D ). (46) 

While one could use this probability to refresh each demon individually, this is 
not very efficient (oc N D steps). We used eq. ( |46D to directly choose new random 
M.D- Depending on the old value of J^Ad, we then either added or subtracted 
magnetization from randomly selected demons, until the right demon magnetization 
was reached. This whole process takes, on the average, only oc \/Nr> steps, and yields 
a new demon magnetization independent of the old one. In our simulations we used 
four times as many demons as spins (iV/5 = 4A^). 

• In order to ensure the canonical energy distribution the spin system is con- 
nected to a heat bath. To make the multispin coding easier, we used a second 
set of demons, this time carrying energy. Individual demon energies vary in units 
of 4: 0,4,8,.... Before each update sweep through the spin system, every spin 
is connected to a magnetic demon and an energy demon; the demons are chosen 
with random order to ensure fast mixing. A spin flip is accepted if and only if the 
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demons can absorb the change in energy and magnetization. The energy demons 
are periodically refreshed with a heat bath. Note that one could also perform a 
normal Metropolis or heat bath update without the energy demons; the update is 
accepted/rejected with the magnetic demon. 

Because the spin update is very fast, we interleaved 5 sweeps through the lattice 
for one demon update. For the largest lattice (128 x 64), an individual spin update 
took 16 ns, whereas the total time divided by the number of spin updates was 28 ns 
on a Cray X-MP. 



Appendix C: Observables with Improved Overlap 

In a 'conventional' correlation length measurement one studies the exponential decay 
of correlators < A t A t + T >, where A t is an observable that depends on the config- 
uration of a single time slice t. On an infinitely long lattice, the correlator can be 
written as 

(0\AT T A\0) 

< A t A t+T >= , - ' . (47) 
+ (0|T T |0) v ; 

Expanding in terms of eigenf unctions of the transfer matrix one obtains 

A|0) = 5>|i>, (48) 



where a, = (z|A|0). Inserting this in eq. (f§7), one obtains 



< A t A t+T >= = £ a 2 eX p(-m J r) , (49) 

A o i 

where the masses are defined by = — ln(^). For a general observable A the a, 
will be nonzero for all i. 

In order to improve this situation one could think of constructing observables 
where only one a, is non-zero. First recall that observables correspond to operators 
which are diagonal in the basis of configuration states \<p). This means that 

(MW) = (<p\A\<P) cW = ^ <W ■ (50) 

The condition that should be the only non-vanishing vacuum overlap is equivalent 
to A\0) = or 

The problem with this equation is that we do not know the wave functions (4>\i) and 
(0| 0) exactly. 

However, we know the eigenfunctions \i)) of the effective transfer matrix. Let us 
embed them into the Hilbert space of the full model. We define 

K)=EI M )(( M K))- (52) 

M 
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Note that (compare section |3|) this projection leaves the vacuum invariant, i.e. 
(i'\0) = ((i\0))- As improved observable we now consider the ratio 



A " ~ (010) ■ (53) 

One can easily convince oneself that A'± can be expressed in terms of the effective 
wavef unctions as 

_ mm 

^-((M(0)|O))- (54) 



The relation of our approach with that of |T(| becomes apparent when one identifies: 

1) Our 5<$, u m with Kronfeld's in his eq. (4.1), where the label M corresponds 
to {€). Note that our correlation matrix is already diagonal for distance 0. 

2) Our A'^ corresponds to Kronfeld's z$ defined in eq. (4.6). 

3) The diagonalization of the effective transfer matrix has its counterpart in the 
variation of C n ^ r (t). 
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Figure Captions 



Figure 1 Overlaps q = Y^m{M\i) 2 for the three lowest states of the 2-dimensional 
Ising model with N a = 6, 8 and 10. 

Figure 2 Eigenstates of the effective transfer matrix T^g for I = 5, 10 and 20. N a 
is 64, and (3 = (3 C . The full line is the ground-state, the dotted line is [0, a)) , the 
dashed line shows |l,s)), and the dash-dotted line gives 1 1, a)). 

Figure 3 Comparison of the convergence of the 'effective correlation lengths' (q/, 
'improved' and 'standard'. 

Figure 4 Comparison of the convergence of the 'effective correlation lengths' 
'improved' and 'standard'. 

Figure 5 The tunnelling correlation length £ ,a> measured from the N a = 64 lattices. 

Figure 6 Eigenstates of the effective transfer matrix T^g for I = 5, 10 and 20. N a 
is 64, and (3 = 0.47. The full line is the ground-state, the dotted line is [0, a)) , the 
dashed line shows |1, s)), and the dash-dotted line gives |1, a)). 
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Table 1: as obtained from the exact effective transfer matrix for N a = 4, using the 
iV r = oo rule 



p 


I 






£l,a 


?2,s 


0.30 


1 


1.50263 


0.51838 


0.28378 


0.21878 


0.30 


2 


1.50263 


0.51880 


0.28378 


0.23855 


0.30 


4 


1.50263 


0.51908 


0.28378 


0.26000 


0.30 


8 


1.50263 


0.51924 


0.28379 


0.28450 


0.30 


4,8 


1.50263 


0.51940 


0.28379 


0.30900 


0.30 


e 


1.50263 


0.51940 


0.28378 


0.28868 


0.45 


1 


5.41722 


0.66675 


0.35470 


0.27231 


0.45 


2 


5.41722 


0.66726 


0.35470 


0.28238 


0.45 


4 


5.41722 


0.66756 


0.35470 


0.29141 


0.45 


8 


5.41722 


0.66771 


0.35470 


0.29822 


0.45 


4,8 


5.41722 


0.66786 


0.35470 


0.30503 


0.45 


e 


5.41722 


0.66786 


0.35470 


0.30243 


0.60 


1 


26.11653 


0.53669 


0.34962 


0.28156 


0.60 


2 


26.11653 


0.53705 


0.34962 


0.28537 


0.60 


4 


26.11653 


0.53726 


0.34962 


0.28802 


0.60 


8 


26.11653 


0.53736 


0.34962 


0.28907 


0.60 


4,8 


26.11653 


0.53746 


0.34962 


0.29012 


0.60 


e 


26.11653 


0.53747 


0.34962 


0.29087 
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Table 2: Correlation lengths £^ as obtained from the exact effective transfer matrix 
for N a = 8, using the N T = oo rule 



p 


I 






£l,a 


6,s 


0.30 


1 


1.57471 


0.67038 


0.38527 


0.26235 


0.30 


2 


1.57479 


0.67318 


0.38709 


0.27564 


0.30 


4 


1.57484 


0.67574 


0.38861 


0.31248 


0.30 


8 


1.57486 


0.67740 


0.38958 


0.35499 


0.30 


4,8 


1.57488 


0.67906 


0.39055 


0.39750 


0.30 


e 


1.57488 


0.67913 


0.39058 


0.41348 


0.45 


1 


12.25453 


1.25922 


0.62386 


0.40659 


0.45 


2 


12.25990 


1.26925 


0.62784 


0.41963 


0.45 


4 


12.26276 


1.27622 


0.63047 


0.43690 


0.45 


8 


12.26420 


1.28000 


0.63191 


0.45234 


0.45 


4,8 


12.26564 


1.28378 


0.63335 


0.46796 


0.45 


e 


12.26565 


1.28381 


0.63338 


0.47114 


0.60 


1 


418.95213 


0.70382 


0.52749 


0.39890 


0.60 


2 


419.01927 


0.71339 


0.53172 


0.40984 


0.60 


4 


419.05376 


0.71861 


0.53417 


0.41564 


0.60 


8 


419.07105 


0.72131 


0.53546 


0.41864 


0.60 


4,8 


419.08834 


0.72401 


0.53675 


0.42164 


0.60 


e 


419.08835 


0.72404 


0.53678 


0.42169 
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Tabic 3: Correlation lengths £^ as obtained from the exact effective transfer matrix 
for N a = 4, using the symmetric periodic lattice rule 





I 


£o,a 




£l,a 


6,s 


0.30 


1 


1.50167 


0.51944 


0.28375 


0.22264 


0.30 


2 


1.50263 


0.51880 


0.28378 


0.23901 


0.30 


4 


1.50263 


0.51908 


0.28378 


0.26007 


0.30 


8 


1.50263 


0.51924 


0.28378 


0.28203 


0.30 


4,8 


1.50263 


0.51940 


0.28378 


0.30399 


0.30 


e 


1.50263 


0.51940 


0.28378 


0.28868 


0.45 


1 


5.41132 


0.66909 


0.35467 


0.27595 


0.45 


2 


5.41713 


0.66737 


0.35470 


0.28348 


0.45 


4 


5.41722 


0.66756 


0.35470 


0.29205 


0.45 


8 


5.41722 


0.66771 


0.35470 


0.29564 


0.45 


4,8 


5.41722 


0.66786 


0.35470 


0.29923 


0.45 


e 


5.41722 


0.66786 


0.35470 


0.30243 


0.60 


1 


26.10328 


0.53833 


0.34962 


0.28395 


0.60 


2 


26.11621 


0.53726 


0.34962 


0.28641 


0.60 


4 


26.11653 


0.53726 


0.34962 


0.28834 


0.60 


8 


26.11653 


0.53736 


0.34962 


0.29040 


0.60 


4,8 


26.11653 


0.53746 


0.34962 


0.29246 


0.60 


e 


26.11653 


0.53747 


0.34962 


0.29087 
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Table 4: Correlation lengths £® as obtained from the exact effective transfer matrix 
for N a = 8, using the symmetric periodic lattice rule 



p 


I 




Cm 


£l,a 


£,2,s 


0.30 


1 


1.56250 


0.67186 


0.26386 


0.38557 


0.30 


2 


1.57475 


0.67310 


0.27658 


0.38685 


0.30 


4 


1.57484 


0.67573 


0.31295 


0.38815 


0.30 


8 


1.57486 


0.67739 


0.35092 


0.38818 


0.30 


A O 

4,8 


1.57488 


0.67905 


r\ o o o on 

0.38889 


0.38821 


0.30 


e 


1.57488 


0.67913 


0.39058 


0.41348 


0.45 


1 


11.68811 


1.29494 


0.41267 


0.63080 


0.45 


2 


12.22444 


1.27367 


0.42339 


0.62794 


0.45 


4 


12.25978 


1.27677 


0.44527 


0.62700 


0.45 


8 


12.26354 


1.27856 


0.52313 


0.65118 


0.45 


4,8 


12.26730 


1.28035 


0.60099 


0.67536 


0.45 


e 


12.26565 


1.28381 


0.63338 


0.47114 


0.60 


1 


404.80434 


0.73313 


0.53976 


0.40926 


0.60 


2 


417.52452 


0.72285 


0.53477 


0.41495 


0.60 


4 


418.82046 


0.72058 


0.53259 


0.42259 


0.60 


8 


418.96946 


0.71637 


0.58641 


0.53098 


0.60 


4,8 


419.11846 


0.71216 


0.64023 


0.63937 


0.60 


e 


419.08835 


0.72404 


0.53678 


0.42169 



Table 5: Estimates £® for N a = 16, iV T = 512. The effective transfer matrix was 
determined by Monte Carlo, using the N T — oo rule. The coupling is (3 = f3 c 



I 






Cl,a 




1 


20.10(8) 


2.402(4) 


1.147(2) 


0.704(1) 


2 


20.17(7) 


2.451(4) 


1.171(3) 


0.726(2) 


3 


20.22(7) 


2.480(6) 


1.181(5) 


0.756(8) 


4 


20.25(7) 


2.499(7) 


1.193(7) 


0.76(1) 


5 


20.26(7) 


2.509(10) 


1.190(15) 




6 


20.27(7) 


2.511(12) 






7 


20.28(7) 


2.518(15) 






8 


20.29(8) 


2.515(20) 






e 


20.339 


2.555 


1.213 


0.8730 


1,2 


20.25(7) 


2.502(6) 


1.195(5) 


0.749(6) 


2,4 


20.32(8) 


2.55(1) 


1.215(15) 




3,6 


20.32(9) 


2.54(2) 






4,8 


20.34(9) 
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Table 6: Estimates £ w for N a = 32, iV r = 1024. The effective transfer matrix was 
determined by Monte Carlo, using the N T = oo rule. The coupling is (3 = (3 C 



I 


?0,a 




Cl,a 




1 


39.28(12) 


4.44(1) 


2.099(2) 


1.266(2) 


2 


39.86(13) 


4.61(1) 


2.176(3) 


1.316(2) 


3 


40.15(14) 


4.71(1) 


2.225(5) 


1.350(4) 


4 


40.32(15) 


4.77(1) 


2.259(6) 


1.385(5) 


5 


40.45(16) 


4.83(1) 


2.283(8) 


1.414(6) 


6 


40.56(17) 


4.87(1) 


2.31(1) 


1.445(9) 


8 


40.68(18) 


4.92(2) 


2.34(2) 




10 


40.75(18) 


4.96(2) 


2.34(3) 




e 


40.727 


5.097 


2.404 


1.7014 


1,2 


40.45(15) 


4.78(1) 


2.258(6) 


1.371(4) 


2,4 


40.79(17) 


4.96(2) 


2.348(11) 


1.461(9) 


3,6 


40.97(20) 


5.05(2) 


2.395(20) 


1.55(2) 


4,8 


41.05(22) 


5.08(2) 


2.42(4) 




5,10 


41.05(22) 


5.10(3) 
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Table 7: Estimates £ w for N a = 64, N T = 2048. The effective transfer matrix was 
determined by Monte Carlo, using the N T = oo rule. The coupling is (3 = (3 C 



I 






£l,a 




1 


73.62(31) 


8.13(1) 


3.783(5) 


2.271(3) 


2 


76.00(31) 


8.54(2) 


3.985(4) 


2.396(5) 


3 


77.32(32) 


8.82(2) 


4.119(6) 


2.483(4) 


4 


78.15(31) 


9.01(2) 


4.215(7) 


2.547(5) 


5 


78.74(33) 


9.15(2) 


4.291(8) 


2.598(7) 


6 


79.17(34) 


9.27(2) 


4.346(9) 


2.641(9) 


>-7 

7 


79.50(34j 


9.36(2) 


A on/ 1 \ 

4.39(1) 


2.671(1) 


o 


79.77(34) 


n A A /o\ 

9.44(2) 


4.43(1) 


2.70(1) 


y 


79.98(35) 


9.50(2) 


4.46(1) 


2.73(1) 


i n 

10 


80.13(36) 


9.56(2) 


4.49(2) 


2.76(2) 


1 o 

L2 


80.37(36) 


9.64(2j 


A A ( 

4.54(2] 
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on rr /i / o'7\ 

80.54(37) 


9.70(3) 


4.57(2) 




lb 


on ft - ? ( o o \ 

80.67(38) 


9.75(3) 


4.61(3) 




on 

2(J 


on o o / 1 n \ 

80.88(40) 


n o 1 / a \ 

9.81(4) 






e 


81.479 


10.188 


4.797 


3.4014 


1,2 


78.53(32) 


9.00(2) 


4.210(8) 


2.536(5) 


2,4 


80.43(35) 


9.52(3) 


4.47(1) 


2.72(1) 


3,6 


81.12(40) 


9.76(3) 


4.60(2) 


2.82(2) 


4,8 


81.46(43) 


9.91(3) 


4.66(2) 


2.88(2) 


5,10 


81.57(42) 


10.00(3) 


4.72(3) 


2.95(3) 


6,12 


81.59(48) 


10.05(4) 


4.75(3) 




7,14 


81.62(43) 


10.08(4) 


4.77(5) 




8,16 


81.59(43) 


10.09(5) 






10,20 


81.64(46) 


10.08(7) 







Table 8: 'Conventional' correlation length estimates £o,a* an d Ci,/ sd, /3 — /3 C from 
two-slice correlation functions. The lattice is N a = 64 by N T = 2048 



t 


10 


20 


40 


80 








*ett",t 
S0,a 


78.9(5) 


81.4(6) 


82.0(8) 


83.4(14) 








t 


1 


3 


5 


9 


11 


13 


15 


f ett",t 

Sl,s 


6.95(2) 


8.33(2) 


9.05(4) 


9.68(6) 


9.81(7) 


9.95(8) 


9.95(11) 
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Table 9: Statistics stat for the multimagnetical simulations at (3 — 0.47, given in 
units of one million lattice sweeps 





16 


16 


32 


32 


32 


64 


64 


64 




16 


32 


16 


32 


64 


32 


64 


128 


stat 


20 


40 


40 


50 


40 


50 


50 


60 



Table 10: Estimates £® for (3 = 0.47 and N a = 16, obtained from the effective 
transfer matrix with the asymmetric rule 





I 






£l,a 


16 


1 


76.3(7) 


2.041(5) 


1.131(2) 


16 


2 


76.5(8) 


2.104(5) 


1.158(2) 


16 


4 


76.7(8) 


2.163(6) 


1.173(6) 


16 


8 


76.6(8) 


2.21(1) 




32 


1 


78.5(4) 


2.013(4) 


1.126(1) 


32 


2 


78.8(4) 


2.076(4) 


1.155(2) 


32 


4 


78.9(5) 


2.136(5) 


1.179(4) 


32 


8 


79.0(5) 


2.18(1) 




32 


16 


79.1(5) 






e 


oo 


78.159 


2.205 


1.218 



Table 11: Estimates for (3 = 0.47 and N a = 32, obtained from the effective 
transfer matrix with the asymmetric rule 



N T 


I 






£l,a 


16 


1 


558(9) 


3.42(2) 


1.914(6) 


16 


2 


607(10) 


3.53(2) 


2.034(6) 


16 


4 


644(11) 


3.65(2) 


2.172(8) 


16 


8 


656(13) 


3.72(2) 


2.27(1) 


32 


1 


718(9) 


2.68(1) 


1.779(4) 


32 


2 


728(10) 


2.84(1) 


1.879(5) 


32 


4 


730(10) 


3.02(1) 


1.982(5) 


32 


8 


731(10) 


3.18(1) 


2.03(1) 


32 


16 


732(11) 


3.30(2) 


2.15(1) 


64 


1 


730(7) 


2.65(1) 


1.776(3) 


64 


2 


740(7) 


2.80(1) 


1.875(3) 


64 


4 


746(8) 


2.98(1) 


1.980(4) 


64 


8 


751(8) 


3.13(1) 


2.06(1) 


64 


16 


757(9) 






64 


32 


760(10) 






e 


oo 


753.48 


3.311 


2.198 
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Table 12: Estimates £ w for (3 = 0.47 and N a = 64, obtained from the effective 
transfer matrix with the asymmetric rule 



N T 


I 






Cl,a 


32 


2 


25800(900) 


4.00(4) 


2.52(2) 


32 


4 


30600(900) 


4.15(4) 


2.80(2) 


32 


8 


33800(1100) 


4.29(4) 


3.12(2) 


32 


16 


35400(1900) 


4.38(4) 


3.35(3) 


64 


4 


38700(2100) 


3.15(6) 


2.60(3) 


64 


8 


38600(2000) 


3.46(5) 


2.84(4) 


64 


16 


39500(2000) 


3.72(5) 




64 


32 


41000(2200) 






128 


16 


41500(1400) 






128 


32 


42800(1500) 






128 


64 


43500(1600) 






e 


oo 


44014.4 


4.002 


3.311 
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